function [varargout] = synfit()

% default values for the table

global CONTROL

sf = getmainselection();


table.F = 0.17; % release probability (constant; no facilitiation)
table.k0 = 1/2.75; % /s, baseline recovery rate from depletion (slow rate)
table.kmax = 1/0.046; % /s, maximal recovery rate from depletion (fast rate)

table.td = 0.015; %  time constant for calcium-dependent recovery
table.kd =  0.09; % affinity of fast recovery proces for calcium sensor

table.ts = 3.5; % decay time constant of glutatme clearance
table.ks = 0.6; % affinity of receptor desensitization for glutatmate
% The large value means no desense occurs (otherwise, ks should be about
% 0.6)

table.kf = 1.5; % affinity of facilitiation process
table.tf = 0.018; % make facilitation VERY slow

table.dD = 0.37; % sets Ca that drives recovery(Ca influx per AP)
% 0.02 yields rate-dep recovery in 100-300 Hz
table.dF = 0.07; % sets Ca that drives facilitation

table.glu = 0.4;

f = fieldnames(table);
p = zeros(length(f), 1);
for i = 1:length(f)
    p(i) = eval(['table.' char(f{i}) ';']);
end;


lb = [0.05,  table.k0/10,  table.kmax/5, ...
    table.td/20,  table.kd/100,  ...
    table.ts/20,  table.ks/20, ...
    table.kf/5,  table.tf/2, ...
    table.dD/5, table.dF/5, 0.001];
ub = [1,     table.k0*10, table.kmax*5, ...
    table.td*20, table.kd*10, ...
    table.ts*10, table.ks*80, ...
    table.kf*20, 1, ...
    table.dD*10, 3, 20];
options = optimset('LargeScale', 'off', 'Algorithm', 'active-set');
[x, fval] = ...
    fmincon(@depfits2, p, [], [], [], [], lb, ub, @synconstraint, options); %#ok<NASGU>
for i = 1:length(x)
    fprintf(1, '%-12s: %.6f\n', char(f{i}), x(i));
    eval(sprintf('fittab.%s=%f;', char(f{i}), x(i)));
end;
XuF(1, fittab);
if(nargout > 0)
    varargout{1} = fittab; % returns the results
end;

% provide evaluation of error sensitivity for each parameter
j = 1;
tb = cell(length(sf), 1);
ib = tb;
trb = tb;
rb = tb;
ntb = tb;
nr = tb;
for i = sf % get the data sets we need.
    tb{j} = 0.001*CONTROL(i).EPSC_Train.stim;
    ib{j} = mean(CONTROL(i).EPSC_Train.amps, 1);
    trb{j} = 0.001*(CONTROL(i).EPSC_Recovery.Tdelays-CONTROL(i).EPSC_Recovery.rdelay);
    rb{j} = mean(CONTROL(i).EPSC_Recovery.Iamps, 2)';
    j = j + 1;
end;
for i = 1:length(tb)
    if(~isempty(trb))
        ntb{i} = 0.001*([tb{i} max(tb{i})+trb{i}]);
        nr{i} = [ib{i} rb{i}];
    else
        ntb{i} = 0.001*([tb{i} tb{i}(end)]);
        nr{i} = [ib{i} ib{i}(end)];
    end;
    nr{i} = nr{i}/nr{i}(1);
end;

% [xo, yo] = XuF(2, table, tb, trb);


test  = 1.2; % 3 x larger and 3 x smaller... 

f = fieldnames(fittab);
for i = 1:length(f) % for each field in the table

    w  = eval(['fittab.' char(f{i}) ';']); %#ok<NASGU>
    be = zeros(3,1);
    for j = 1:3
        tbtest = fittab; % reset to original table
        switch j
            case 1
                eval(sprintf('tbtest.%s=w;', char(f{i})));
            case 2
               eval(sprintf('tbtest.%s=w*test;', char(f{i})));
            case 3
               eval(sprintf('tbtest.%s=w/test;', char(f{i})));
            otherwise
        end;
        [xo, yo] = XuF(2, tbtest, tb, trb);
        err = 0;
        for k = 1:length(xo)
            err = err + sum((yo{k}' - nr{k}).^2);
        end;
        be(j) = err;
    end;
    fracerr = ((be(2)+be(3))/2)/be(1); % fractional error - average of 2x and 0.5x, relative to best value.
    fprintf(1, '%5s = %f,  errp = %f\n', char(f{i}), eval(['fittab.' char(f{i})]), fracerr);
end;

